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ABSTRACT 



Aims. We have developed a time-dependent two-component hydrodynamics code to simulate radiatively-driven stellar winds from hot stars. 
Methods. We use a time-explicit van Leer scheme to solve the hydrodynamic equations of a two-component stellar wind. Dynamical friction 
due to Coulomb collisions between the passive bulk plasma and the line-scattering ions is treated by a time-implicit, semi-analytic method 
using a polynomial fit to the Chandrasekhar function. This gives stable results despite the stiffness of the problem. 

Results. This method was applied to model stars with winds that are both poorly and well-coupled. While for the former case we reproduce 
the mCAK solution, for the latter case our solution leads to wind decoupling. 
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1. Introduction 

Stellar winds from hot stars are often described using the so- 
called CAK theory jCastor. Abbott. & Klei n, 1975), including 
later i mprovements for a fini te stellar disk dPauldrach et al. 



; [19861 lEliend & Abbott', '1986") and the ionisation stratification 
' of the wind (Abbott, 1982). The radiative line force is here pa- 
rameterised by three constants only, k, a, and 6. 

Although qualitative agreement between theory and obser- 
vations was achieved, there were still discrepancies remaining, 
namely the terminal velocity Vpo being too l ow and the mass 



ture of th e wind s i gnifica ntly. As was first discussed by 
fcassinelli & Olson L [19791) and later by ( Pauldrach & PuTsl 
Il990l) . the X-rays lead to photoinisation via Auger processes 
and are necessary to explain the presence of superionised ions 
O VI in observations. 

Lucv & Solomon d 19701) suggested a new hydrodynamic 



1982; Friend & Abbott, 



loss rate M being too high ("Abbot] 

Ll986). Besides this, ca lculations did not reprod uce the ob 



served ionisation ratios dPauldrach & Puis . 119901) . Further de- 
velopment of the theory was focused mainly on the elTect of 
multiline scat tering, magnetic field an d rapidly rotating stars. 
As shown by Friend & Castor ( 1983b . if one includes the ef- 
fect of multiple scatterings in overlapping lines, the terminal 
velocity is higher than the CAK value, while the mass-loss rate 
remains similar to the latter. Furthermore, the inclusion of a 
magnetic field and the assumption of high rotational velocity 
lead to a higher terminal vel ocity and a higher mass-loss rate 
dFriend & MacGregoiilTgsl . 

The above problem with the ionisation structure can be 
resolved by including X-rays. X-rays are a common fea- 
ture among the O stars and can change the ionisation struc- 



instability of the radiative driving force. This so-called line- 
driven instability can possibly explain observed variabilities 
of stellar winds, such as variable X-ray emission, discrete ab- 
sorption components, and the appear ance of broad absorp - 
tion troughs in P Cygni line profiles (Prinia & Howarth , Il986l) . 
The contra dictory results from the simp hfied hnear stab ihty 
analysis of 'MacGregor et al.' ('1979*) and lCarlber^ dl980l) . on 
one hand, and Abbott (1980), on the other, were unified by 



Owocki & Rvbickil dl984ll985l) in the so-called bridging law. 



The first time-dependent, num erical sim ulations of the instabil- 
ity were undertaken by Owocki et al. ([l988), who found that 
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the wind develop s a train of strong reverse shocks. Their model 
was improved by lFeldmeieii dl995h . including the energy equa- 
tion without the approximation of gas isothermality to calculate 
the temperature structure of the ste llar wind. Adopting asim- 
ple formulation for the line force, iDessart & Owockil d2005h 
extended hydrodynamical instability simulations from one to 
two dimensions (2D), to find the lateral coherence length of 
instability-generated shells of dense gas in the wind. 

The approximation of the one-component flow assumed 
in the CAK model is acceptable for sufficiently dense winds 
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from O stars and hot B stars dCastor. Abbott. & Kleiriill976h . 
The radiative line force acts only on ions (termed "metals" in 
the following) that scatter photons in numerous spectral lines. 
These ions share their momentum through Coulomb collisions 
with the passive part of the plasma (protons). The dynamical 
effect of Coulomb collisions on the plasma is describe d well 
byjiynamical friction, which was introduced by Chandr asekhar 
(Il94i for the gra vitationa l force and later applied to electro- 
magnetic forces bv lSpitzer 

For thin winds from B near-main-sequence stars, however, 
dec oupUng and a well-known plasm a runaway effect can oc- 



1992h . At low wind densities. 
Coulomb interactions are weak, and the momentum transfer 
from metals to protons becomes inefficient. As a result, the two 
components decouple at a given radius. From this decoupling 
radius on, metals are strongly accelerated by the radiative force, 

whereas the passive plasma is de celerated by gravi ty. 

It came as a surprise that iKrticka & Kubi3 i2000h ob- 
tained a different result from stationary hydrodynamic cal- 
culations. They found a wind solution where the two com- 
ponents remain coupled, but undergo a sudden jump to a 
slow-acceleration branch at a definite radius. This is analo- 
gous to what is found for shal low (and overloaded) solutions 
( Feldmeier & Shlo sman'. '2000^ The issue has not been set- 
tled yet ever since Owocki & Puis (2002) showed from lin- 
ear stability analysis that the flow should be disrupted by 
ion separation before reaching a slow-acceleration solution. 
[Krticka & Kubat ( 2001b) also found that in a low-density wind, 
frictional heating may be important for the energy balance of 
the whole wind and may increase the temperature of the gas 
significantly. 

An effect similar t o ion decoupling was considered by 
Porter & Skouzal (Il999h . namely shock decoupling. In a low- 
density wind that passes through a shock, the postshock 
gas remains at high temperatures and in a highly ionised 
state because the gas is too rarefied to cool radiatively 
dKrolik & Ravmondl Il985h . In this case the ions responsible 
for line driving are completely stripped and are not accelerated 
further by the radiation field. 

Another inte resting idea in this cont ext is the generation of 
pulsating shells dPorter & Skouzal Il999l) . If the wind decouples 
at a radius where the flow speed is still lower than the escape 
speed, the passive plasma falls back to the star, and the inter- 
action with the outflow leads to pulsating shells. However, this 
result is based only on a one-component model, with an artifi- 
cial turn-off of the radiative acceleration. 

To analyse the possible occurrence of ion decoupling and 
pulsating shells in thin stellar winds, it is necessary to develop a 
suitable time-dependent, two-component hydrodynamics code. 
In the present paper we describe our numerical method to sim- 
ulate these outflows, and give basic results on thick and thin 
winds. 



sive plasma. The forces acting are gravity, dynamical friction, 
gas pressure gradients, and, for only the line-scattering ions, 
the radiative line force. The continuity equations are 



dp^ 1 5(rVpfp) 
dpi 



dr 

1 dir^piVj) 
dr 



and the Euler equations are 



dv^ dvp 1 dpp 

dt dr pp dr 

dvi dv\ I dpi 

dt ' dr pi dr 



= ^eif, 



- ^rad-^eff- 



(la) 

(lb) 

(2a) 
(2b) 



Here, p;, Ui, and pi are the density, velocity, and pressure of the 
metals, and pp,Vp, and pp are the density, velocity, and pressure 
of the passive plasma. The effective gravitational acceleration 
is ^eff = -GMt{l - Te)lr^, with gravitational constant G and 
Eddington factor Fp. The frictional force Rpi between metals 
and passive plasma is described below (see Eq.[8]). The above 
system is closed by the equations of state. 



Pp = ^kT/nipPp, 
Pi = ^JkT/miPi. 



(3a) 
(3b) 



The assumption of wind isothermality becomes questionable in 
the presence of strong frictional heating, and we aim to include 
the energy equation in future work. 

2.1. Radiative acceleration 

We write the radi ative acceleration in the form 
dKrticka&KubatlbOOOl Eq.3) 



Tad 



f'ionft 



fin, 



(4) 



with CAK force multipliers k, a, 6. Here, L, is the luminosity 
of the star, Uth the thermal velocity of ions, and cr ^ 0.33 cm^ 
g ' the Thomson opacity due to scattering on electrons. 

The radiative acceleration due to line scattering only acts 
on metal ions, but the CAK force multipliers k, a, 6 were cal- 
culated for a one-component plasma. We account for this by a 
scaling facto r 77 to the radiative force . The value 77 = 0.0127 
calculated bv iKrticka & Kubat (l2000h for solar metallicities is 
adopted here. 

The finite disk correction factor ff\„ is given by 



(ICastoret aUll975h 



/fin(r) 



(1 -Hcr)"+1 -(1 +0-pi) 

a-(a-+ -pl)(a+ 1) 



(5) 



2. Two-component winds 

We restrict ourselves to a ID spherically symmetric, isother- 
mal, quasineutral, two-component outflow consisting of metals 
that scatter stellar photons in numerous spectral lines and a pas- 



where = ^j\-(Rt |l^ (with stellar radius /?»), and cr is 
given by dCaston,ll974l) : 



r dvi 
Vi dr 



(6) 
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Finall y, f\„„ is a corr ection for the ionisation state of the stellar a staggered mesh (see lLeVeque et allll998l p. 131). We calcu 



wind ( lAbbottlll982h . 



fion(r) 



10-"[cm3]ne 
W{r) 



(7) 



where W(r) = ^(1 - ^f^^(R7^¥) 

is the geometrical dilution 
factor Due to the assumption of quasineutrality, the number 
density of electrons roughly matches the number density of 
the passive plasma. Thus, we use - Hp in /ion. This correc- 
tion factor does not hav e a significant influence on the wind 
dynamics (lAbbottl[l982l) . 



2.2. Friction terms 

The passive plasma and absorbing ions interact via Coulomb 
colhs ions, which are described by a fr ictional force Rpi per vol- 
ume dSpringmann & Pauldrachlll992h . 



(8) 



where rip and rii are the number densities of the passive plasma 
and absorbing ions, respectively, and the frictional coeflicient 
^pi is given by 



^pi 



kjiT h-Up\' 



(9) 



Here, Zjg and Zpg are the ion and passive plasma charges, re- 
spectively. The Coulomb logarithm In A is defined as 



hiA^ln 



2471 1 keT \ 



1.5 



(10) 



with Boltzmann constant ks, total number density n, and 
wind temperature T. The Chandras ekhar function G(xpi) 
( Chandrasekhar , 1943 : Spitzer , 19561) in Ea.(|8]l i s given in 
terms of the error function <l>(xpi) by ( Spitzei , ll956h 



G(jCpi) = 



1>(xpi) exp(-xgi) 

Xpi 



pi 



(11) 



This function depends on jCpj, the ion separation drift speed rel- 
ative to th e passive plasma, scaled to the m ass- weighted ther- 
mal speed jSpringmann & Pauldrachll 19921) 



l^i - Vp\ 



wth VI +MIAp 



(12) 



where Ai and Ap represent the mean atomic mass of ions and 
passive plasma in atomic units. 

3. Method of solution 

To solve the four hydrodynam ic equations we u se a hydro- 
dynamics code developed by iFeldmeieii (119951) as a core. 
It employs a standard Euler scheme and is suited for ID, 
one-component outflows. We extended this code to a two- 
component version. Equations ([TJ and (|2]i are discretised using 
an operator-splitting, time-explicit, finite difference method on 



late advection fluxes using van Leer's monotonic interpolation 

I 1 1 — — ! 

(see'van Leer, 1982). 



Krticka & Kubatj( 2000l) discussed the circumstances under 
which the radiative line acceleration is balanced by dynamical 
friction and not by inertia. In this case, a decrease in dynamical 
friction leads (counter-intuitively) to a decrease in the radiative 
acceleration of the gas. We accelerate the gas by the sum of dy- 
namical friction and radiative force, i.e. avoid operator-splitting 
of these two force terms, in order to achieve a stable numerical 
scheme. Similarly, for a barometric density stratification, one 
has to apply the sum of the thermal pressure force and gravity 
at once to avoid an unstable scheme. 

As the time step, we use the minimum of the separate 
Courant time steps for the two flow components. The wind is 
characterised by two parameters, the mass-loss rate and termi- 
nal speed. The latter scales with the escape speed, and the for- 



mer is given roughly by the CAK relation (iCastor et al. 
Eq.46) 



1975 



McAK 



47rGM 



a{l — a) 



(^r,)i/«(i - r,) 



-(\-a)la 



(13) 



To calculate the radiative force we use tabulated values of the 



CAK parameters from lAbbottI (|1982|) . 

As initial conditions for the case of a well-coupled wind 
we use, for the velocity of the both components in the subsonic 
part, v{r) = 0. la exp {Hr), where H is the scale height and a the 
isothermal sound speed. In the supersonic part, we use v{r) = 
a + 2Qa r. Initial values for the density of both components are 
calculated from Eqs. ([1]) and ( fT3] ). 

For the case of a low-density wind (where poor coupling 
between the components may be expected), we artificially in- 
crease the frictional coupling by increasing the average ion 
charge qi and run the simulation (with the above initial con- 
ditions) until a converged CAK one-component solution is 
achieved. This solution then serves as initial conditions for both 
components in a subsequent simulation with a realistic average 
ion charge. 

Our bound ary conditions are set according to the theory of 
characteristics (lAndersonL 1 19951 pp. 303-307). For absorbing 
ions and passive plasma, we keep the densities pp and pi at the 
inner boundary fixed to their initial values and set, at each time 
step, the momentum densities piVi, p^Vp at the inner boundary 
to their value at the first interior mesh point (zero-order extrap- 
olation). At the outer boundary, we extrapolate the mass and 
momentum densities of both components from the last interior 
mesh point. 

3.1. Discretisation of friction 

The central issue of the present paper is the inclusion of the 
frictional term (O. In the following we develop a simple, com- 
putationally effective, and accurate method for determining of 
the velocity difference Xpi at every time step, based on the fric- 
tional term. 

The equations for the momentum change of the passive 
plasma and the absorbing ions due to the frictional term alone 
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Fig. 1. The Chandrasekhar function; exact function G(Xpi) ac- 
cording to ( fTTI ) (dashed line), and approximation GA(JCpi) (solid 
line). 



are, from (O 

dVrJ _ CG(Xpi) 

dt pp 

dVr,p _ CG(Xpi) 

dt Pi 
where the constant C is given by 

C 



(14) 
(15) 



AiApmj 



2 ■ 



Transforming the above differential equations into time- 
explicit difference equations results in a system with high stiff- 
ness, since we attempt to evolve a stellar wind on a relatively 
slow flow time scale with respect to which the faster, frictional 
processes maintain equilibrium almost instantaneously. As a 
consequence, the numerical solution fails due to large oscilla- 
tions at the wind base. Therefore we set up the following semi- 
implicit scheme for the frictional terms. Subtracting Eq. ( fTSl l 
from Eq. ( fT4l i. we obtain a differential equation for the drift ve- 
locity. 



dxpi 



(pi + pp)C 



:G(Xpi). 



(17) 



fth VI +A-JAp 

In Eqs. (fT4l i and ( fTsT i we now replace the partial time derivatives 
by a total derivative following the idea of operator-splitting, i.e. 
,friction is now treated without reference to any "forthcoming" 
hydrodynamic processes. The djdt is not to be understood as a 
Lagrange derivative. Upon formal integration. 



pi ' d.xpi ^ _ r 



to V(a Vl + A 



(18) 



where A denotes the ratio A JA^, and Xpj o, resp. jCpj i, are the 
drift velocities at times fo, resp. to + Af. We want to obtain an 
analytical expression for the left hand side of Eq. ( fTSl l. To this 
end we use a three-interval fit to the Chandrasekhar function. 



For low and high drift velociti es we use the same approxima- 
tions as lOwocki & Pulsl (|2002|) . 



and 



Gi(xpi) 



G^ixpi) 



7x ■ 
3V^ 



forxpi < xi 



K 



lx\ 
pi 



2 for.Xpi>X2. 



(19) 



(20) 



For values between x\ and x^, we approximate the 
Chandrasekhar function using a quadratic function. 



Gi{xp^ » ajXpi + b2Xpi + C2, 



(21) 



where the parameters a2,b2, C2 are evaluated from the condi- 
tions Gi(xi) = G2(xi), G2(xinax) - 0.214, and dGi(xi)/dx - 
dG2(xi)/dx. Continuity of the derivative at X2 is not required, 
because we have only a second-order polynomial approxima- 
tion. The points xi , X2, and the scaling factor K were chosen 
to achieve the best fit to the Chandrasekhar function (fTTl i. We 
find, for xi = 0.1 and X2 - 1.5, the values 02 - -0.2341, 
b2 = 0.4532, C2 = -0.0053 and K = 0.74. This approximation 
to the Chandrasekhar function is termed Ga in the following 
and is shown in Fig[T] 

For a well-coupled wind the difference between G and Ga 
is of minor importance, and both functions lead to essentially 
the same wind solution. To show this, we use two different ap- 
proximations Ga, one being a global overestimate, the other a 
global underestimate of G. In both cases, the code converged 
to almost the same steady CAK solution. 

Using the approximation Ga, we obtain the velocity differ- 

(16) 

ence Xpi as follows. At the dynamically most important part of 
the Chandrasekhar function, i.e. its maximum (covered by G2), 
integration of ( fTST i leads to 



In 



2a2JCpi +b2-E 



2a2Xpi + b2 + E 



= -CArS, 



(22) 



.Vpi,0 



where H - ^bl, - 4fl2C2- With the help of the substitutions 



- 2a2Xpifl -H ^2 - S, 
q+ - 2a2Xpio + b2 + S, 



(23) 
(24) 



we can simplify the expression (l22l l and as final form of the 
velocity difference obtain 



1 

2a2 



1 -H 2-exp(-CAfH) 
1 - 2^exp(-CAfH)'^ 



-H — bo 



(25) 



For low drift velocities and using the approximation Gi(xpi), 
integration of ( fTSl ) gives 



Xpi,i = Xpi,oexp 



3 



CAf 



(26) 



and finally, for high drift velocities and using the approxima- 
tion Gi(xpi), 



^pi,o ■ 



^KCAt. 



(27) 
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The expressions ( |25|26|27] i are used in our code to calculate 
the change in drift velocity Xpi due to friction during a hydro- 
dynamic time step. This results in a highly improved stabil- 
ity behaviour compared to direct, time-explicit differencing of 
the friction terms and, for the first time, allows time-dependent 
simulations of multi-component winds. 

Finally, we note that a co mmon ly used approximation to 
the Chandrasekhar function, (lKarlickvetalll2000h 

which covers the whole interval of drift speeds, results in inte- 
grals we were not able to carry out analytically. 

3.2. Time requirements 

Since the two-component code is essentially a duplication of 
the one-component code, and since the frictional force terms 
are solved by analytic equations, i.e. at very small computa- 
tional cost, the cpu time for a simulation is comparable to the 
one-component case, typically a few hours (for 0{lO^) mesh 
points) on a dedicated workstation. However, we chose a very 
small Courant number of 0.05 (instead of typically 0.5) to bring 
the inviscid Courant time step somewhat closer to the frictional 
time step, in order to avoid changes too large in friction terms 
during subsequent time steps, which could trigger instability. 

4. Results of calculations 

In the present paper , we co nsidered essentially the same stars 
as Krticka & Kubal ( 2000l) . More specifically, we apply our 



method to the wind where a well-coupled solution is predicted 
(we refer to this case as BO) and to the wind with possible de- 
coupling (referred to as B5). The corresponding stellar param- 
eters are giv en in Table [T] Details about these parameters are 
discussed by lKrticka & Kubat (2000). 

To determine whether our model winds are well or poorly 
coupled, we considered the parameter Tl, which is the ratio of 
the radiative force p er mass acting on ions, g' , , to th e gravita- 
tional acceleration g ISpringmann & PauldrachTd 19921) . 



'rad 



'rad 



(29) 



Pi+Pv S 8 
and the parameter Fb, which is the ratio of the CAK radia- 
tive force g'^^J^ to gravity, for the case that dynamical friction 
reaches a maximum. 



„CAK\ 



l.i„i=0.%8 



(30) 
8 

Decoupling should occur, if 

FB<rL. (31) 
^From Fig|2] decoupling is expected in model B5 for r > 1 .6Rt 



Analy tically, is determined by ISpringmann & Pauldrach 
(Il992h . 

l/(3/?-l)\-' 



1 - 



/ McakJ] 



where /3 is the beta-law parameter (see, e.g.. |Porter & Skouza . 
1999h and T] has been introduced in (|4|. 



4.1. Well coupled stellar wind 



For this case we used the BO model star from lKrticka & Kubat 
(12000) . Our hydrodynamic simulation evolved to a steady 
mCAK solution with the correct terminal speed v^o- Figure [3] 
shows the velocity law and the densi ty stratification of the 
wind. For comparison with , Krticka & Kubat (2000) and the 
CAK predictions, we summarise the results of our calculation 
in Table |2] As was expected, for the case of a well-coupled 
wind, the velocities of both absorbing ions t;„^, i and passive 
plasma Uoo,p are roughly the same. 

The agreement between our results and those obtained with 
the one-component version of the code is very good, see Fig. [3] 
There is also good agree ment with the stationary calculation by 
Krticka & Kubi^ (I2OOOI) . 

To test the sensitivity of the results to our approximation of 
the Chandrasekhar function, we used different sets of param- 
eters a2,b2,C2 and of points x\,X2 without the condition of a 
best fit. For 02 = -0.176, ^2 = 0.404, C2 = -0.001, xi = 0.15, 
and X2 = 1.85, fit Ga lies everywhere above G, with maxi- 
mum deviation ^ 30%, whereas for 02 - -0.224, b2 - 0.438, 
C2 = -0.0033, xi = 0.10, and X2 = 1.10, the function Ga 
lies everywhere below G, with maximum deviation « 50%. It 
turns out that, for the case of a well-coupled wind, the detailed 
form of the approximation of the Chandrasekhar function is not 
very important, since the steady-state wind solutions are always 
fairly similar. On the other hand, for a poorly coupled wind, a 
correct approximation is mandatory, since the frictional force 
determines the point where the wind starts to decouple. 

4.2. Low-density wind with decoupling 

As a second model, we considered the B5 star from 
Krticka & Kubat (2000). Their two-component model shows a 
solution with lower acceleration compared to the normal CAK 
solution, but did not lead to decoupling. R esults of our cal- 
culations are shown in Fig. H] Compared to iKrticka & Kubat 
(2OOO), we changed the average ion charge to a slightly higher 
value qi = 2.5 instead of qi = 2.0, as before. The aim of this 
was to increase dynamical friction to prevent the appearance of 
pulsating shells. Namely, this higher value of qi shifts the de- 
coupling radius to a location where the escape velocity is lower 
than the local speed of the passive plasma and, consequently, 
matter is no longer gravitationally bound to the star. The more 
subtle case when the matter is still gravitationally bound to the 
star will be considered in a forthcoming paper. 



Table 1. Parameters of wind models 



(32) 



Star 


M 




R, 


a 


k 


S 






[Mo] 


[K] 


[Ro] 










BO 


90.0 


28500 


37.0 


0.590 


0.170 


0.09 


3.0 


B5 


4.36 


15500 


3.01 


0.511 


0.235 


0.12 


2.5 
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Fig. 2. Left panel: Plot of Fb and Fl as function of radius for the model B5. The curves cross at r ~ 1.6/?*. Right panel: 
Dependence of the decoupling radius rd on mass-loss rate M for the model B5. If the wind is denser, the decoupling radius is 
located farther away from the star 
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Fig. 3. Left panel: Final wind-velocity laws for model BO after 55 flow time units Rt/voo- Absorbing ions are marked by ( — 

), passive plasma by ( ), and the result of the one-component code by (■ ■ ■ ). The curves for ions and passive plasma are 

indistinguishable. Right panel: Density for the same model BO. 



Table 2. Calculated values of Voo,p and y„^4 in comparison with 
u'^^ from the one-component model, and derived values of 
the mass-loss rate for absorbing ions, Mi, and for the passive 
plasma, Mp. 



Star 








Mi 


Mp 




[km.s-1] 


[km.s-1] 


[km.s-1] 


[Me p. y.] 


[Mo p. y.] 


BO 


1600 


1600 


1600 


6.6 10-** 


4.4 10-^ 


B5 


200 


36000 


800 


2.3 10-'" 


1.510-'^ 



Note L The value of the passive plasma mass-loss rate roughly cor- 
responds to the CAK value. 



We find that metal ions decouple from the passive plasma 
and start to accelerate steeply at the decoupling radius, whereas 
the passive plasma starts to decelerate at this location. The de- 
coupling point from the simulation agrees roughly with the pre- 



diction from expression ( |32] |. rd » 1.6. (To derive this value, we 
use the CAK value for v^, from Table |2] and (3 - 0.8, which is 
a good estimate for our model.) The present decoupling con- 
tradicts the results of Krticka & Kubat (2000), who obtained a 
shallow, coupled solution. 

The quite un expected res ult of Krtick a & Kubat (l2000l) 
was analysed by lOwocki & Puis (2002) . who performed a 
linear stability analysis of the time-dependent hydrodynamic 
equations to derive perturbation growth rates and propaga- 
tion speeds in a two-component stellar wind. These authors 
found that the i on de coupling instability persists for the 
Krticka & Kubat (2000*) shallow-wind solution, for long wave- 
length perturbations with h igh temporal growth rate ^ lO'^y/i?,. 



It therefore seems that the Krticka & Kubat ( 2000l) solution is 
not a physically valid solution. 



Owocki & Puis (2002) also show that for the solution ob- 



tained by Krticka & Kubat (2000l) . the ion-decoupling instabil- 
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ity has a modest spatial growth rate, mainly due to the fast 
speed of perturbation propagation, which only allows a slight 
amplification during the time needed to converge to the steady- 
state solution. W e may speculate that this is related to how 
Krticka & KubatI (I2OO0I) indeed found convergence of the it- 
eration scheme applied to solve the steady-state equations to a 
shallow, one-component solution. 

As a test, we increased the value of the ion charge to a 
rather unrealistic value of q-, = 5, which artificially increases 
the Coulomb coupling. For this case we obtained a stable, one- 
component flow solution according to CAK. 

After decoupling, the ion velocity gradient is so large that 
the associated spectral lines from ions should be optically thin. 
Because of this dramatic reduction of the Sobolev optical depth 
in ions, the absorption of this highly accelerated material be- 
comes weak and the signature of this material in the spectrum 
is also weak. We must mention that the CAK radiative acceler- 
ation given by Eq. (|4| is overestimated. Babel ( 199 6^ showed 
that inclusion of the shadowing effect by photospheric lines to 
radiative acceleration calculation lowers its value significantly 
for B stars with a low-density wind compared with the CAK 
model. Thus, the high ion speeds from the above figure would 
not be observed, an d the observed wind speed is inst ead that of 
the passive plasma dSpringmann & Pauldrachill992h . 



5. Summary 

We develop a method for time-dependent hydrodynamic simu- 
lations of multi-component stellar winds. To avoid resolving 
the prohibitively short friction time scale (causing very stiff 
system of equations) in our time-explicit scheme, we used a 
three-interval fit to the Chandrasekhar function and solved the 
friction terms (after operator-splitting) analytically. The stabil- 
ity and accuracy of our method is demonstrated for a BO model 
star with a well-coupled wind, where our flow solution evolves 
to the well-known steady mCAK solution. 

For a B5 model star with a low-density wind, we find that 
ion decoupling from the passive wind plasma occurs at a defi- 
nite wind radius, instead of the transition to a shallow, coupled 
wind solutio n with small acceleratio n of both components, as 
predicted bv lKrticka & KubatI ( I2OO0I) . 

In the future, we will apply the code to a larger sample of 
stars and study whether p ulsating shells may origi nate in multi- 
component stellar winds ( Porter & SkouzaL 19991) . 
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Fig. 4. Left panel: Final wind-velocity laws for model B5 after 55 flow time units R,/Voo- Absorbing ions are marked by (• • • ), 
passive plasma by ( — ), and the result of the one-component code by ( — ). Right panel: Wind density for the same model. 



